Niniejszy raport przedstawia szczegółówą analizę właściwości elektrod grafenowych na podstawie zbioru danych. Dane zawierają 925 próbek, z czego znaczna jej część jest niekompletna. Do ich uzupełnienia zastosowano metodę k-nearest neighbors. Następnie przeprowadzono analizę statystyczną, sprawdzono rozkłady wartości atrybutów oraz, na podstawie macierzy korelacji, znaleziono kilka ciekawych zależności.
Jako rezultat przeprowadzonej analizy zaobserwowano:
1. Udział węgla ma istotny wpływ na właściwości
materiału. Zwiększenie zawartości węgla w materiale wysoce
wpływa na zmniejszoną pojemność materiału. W zamian zyskujemy
zwiększenie okna stabilności.
2. Pojemność rośnie wraz z udziałem tlenu. Niski udział
tlenu (<25%) praktycznie gwarantuje niską pojemność materiału, ale
zapewnia wysokie okno stabilności. Udziął azotu ma natomiast marginalny
wpływ na pojemność.
3. Pojemność jest negatywnie skorelowana z oknem
stabilności. Trudno jest znaleźć materiał który ma jednocześnie
wysoką stabilność i ma dużą pojemność.
W części modelowania podjęto próbę przewidzenia pojemności materiału
na podstawie pozostałych informacji ze zbioru. Skorzystano z dwóch
modeli uczenia maszynowego:
1. Model regresji liniowej, który okazał się
niewystarczający, szczególnie dla próbek o wysokiej pojemności,
2. Model Random Forest, który poradził sobie
zdecydowanie lepiej, dla niższych wartości pojemności trafia praktycznie
w punkt, a dla wyższych tylko nieznacznie je niedoszacowuje.
Na zakończenie zastosowano analizę SHAP dla modelu Random Forest pokazując przykładowe rozumowanie algorytmu.
Do przygotowania raportu skorzystano z następujących bibliotek przedstawionych na zajęciach:
library(ggplot2)
library(dplyr)
library(tidyr)
library(ggsci)
library(VIM)
library(tibble)
library(knitr)
library(plotly)
library(randomForest)
library(caret)
library(shapper)
library("DALEX2")
Dla zapewnienia powtarzalności wyników skorzystano z funkcji set.seed()
set.seed(123)
Zbiór danych dotyczy właściwości elektrod grafenowych i zawiera szczegółowe informacje na temat ich budowy, warunków testowania oraz uzyskanej wydajności.
Plik jest wczytywany przy założeniu, że jest umieszczony w tym samym katalogu co plik z tym raportem.
df <- read.csv('data.csv')
df <- as_tibble(df)
Wszystkie kolumny zostały przemianowane na jej skrócone odpowiedniki, aby ułatwić przetwarzanie danych. Usunięto kolumnę ref, ponieważ jest to jedynie identyfikator którego nie będziemy używać przy analizie danych. Dodatkowo, usuwamy pw_limits, z uwagi na to, że powiela ona informację z pw_low i pw_up
df <- df %>%
rename(
ref = `Ref.`,
pw_limits = `Limits.of.Potential.Window..V.`,
pw_low = `Lower.Limit.of.Potential.Window..V.`,
pw_up = `Upper.Limit.of.Potential.Window..V.`,
pw = `Potential.Window..V.`,
current_density = `Current.Density..A.g.`,
capacitance = `Capacitance..F.g.`,
ssa = `Specific.Surface.Area..m.2.g.`,
rct = `Charge.Transfer.Resistance..Rct...ohm.`,
rs = `Equivalent.Series.Resistance..Rs...ohm.`,
electrode_conf = `Electrode.Configuration`,
pore_size = `Pore.Size..nm.`,
pore_volume = `Pore.Volume..cm.3.g.`,
id_ig_ratio = `Ratio.of.ID.IG`,
n_at = `N.at.`,
c_at = `C.at.`,
o_at = `O.at.`,
electrolyte_formula = `Electrolyte.Chemical.Formula`,
electrolyte_cond = `Electrolyte.Ionic.Conductivity`,
electrolyte_conc = `Electrolyte.Concentration..M.`,
cell_conf = `Cell.Configuration..three.two.electrode.system.`
) %>%
select(-ref, -pw_limits)
head(df)
## # A tibble: 6 × 19
## pw_low pw_up pw current_density capacitance ssa rct rs
## <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 0 0.8 0.8 1 680 186. NA 7.7
## 2 0 1 1 1 367 537 6.1 1.95
## 3 0 1 1 2 338 537 6.1 1.95
## 4 0 1 1 5 283 537 6.1 1.95
## 5 0 1 1 10 246 537 6.1 1.95
## 6 0 0.5 0.5 1 872 168. NA 0.8
## # ℹ 11 more variables: electrode_conf <chr>, pore_size <dbl>,
## # pore_volume <dbl>, id_ig_ratio <dbl>, n_at <dbl>, c_at <dbl>, o_at <dbl>,
## # electrolyte_formula <chr>, electrolyte_cond <int>, electrolyte_conc <dbl>,
## # cell_conf <chr>
Poniżej przedstawiona została ilość brakujących danych dla każdej kolumny.
count(df)
## # A tibble: 1 × 1
## n
## <int>
## 1 925
df %>%
summarise(across(everything(), ~ sum(is.na(.)))) %>%
as.data.frame()
## pw_low pw_up pw current_density capacitance ssa rct rs electrode_conf
## 1 4 4 5 16 17 572 786 772 0
## pore_size pore_volume id_ig_ratio n_at c_at o_at electrolyte_formula
## 1 769 729 596 690 699 703 0
## electrolyte_cond electrolyte_conc cell_conf
## 1 99 62 0
Dodatkowo, później przy analizie danych zauważono, że niektóre wartości kolumn tekstowych są puste, lecz w formie pustej treści a nie jako “NA”. Takie przypadki również podmienimy na NA.
df <- df %>%
mutate(across(
.col = where(is.character),
.fns = ~na_if(.x, '')
))
Skala brakujących danych jest poważna. W kilku kolumnach zdecydowana większość obserwacji ma wartości puste. Jako rozwiązanie zaproponowano zastosowanie algorytmu K-nearest-neighbors, aby uzupełnić kolumny z jak najmniejszym możliwym wpływem na statystyki kolumny jak np. mediana wartości.
Dodatkowo, dla kilku wybranych metryk dodamy nową kolumnę zapisującą informację czy dana wartość była początkowo pusta. Zrobiono to, ponieważ bywają sytuację w których brak informacji również jest wartościową informacją (przykładem - stosunek odmów odpowiedzi w ankietach wyborczych jest różny dla sympatyków różnych partii).
df <- df %>%
mutate(rct_na = is.na(rct), rs_na = is.na(rs), n_at_na = is.na(n_at))
df <- kNN(df, k=5, imp_var = FALSE)
Dane po przetworzeniu brakujących danych zawiewrają następujące dane:
Zbiór zawiera 925 próbek.
head(df)
## pw_low pw_up pw current_density capacitance ssa rct rs
## 1 0 0.8 0.8 1 680 186.3 18.5 7.70
## 2 0 1.0 1.0 1 367 537.0 6.1 1.95
## 3 0 1.0 1.0 2 338 537.0 6.1 1.95
## 4 0 1.0 1.0 5 283 537.0 6.1 1.95
## 5 0 1.0 1.0 10 246 537.0 6.1 1.95
## 6 0 0.5 0.5 1 872 168.2 0.6 0.80
## electrode_conf pore_size pore_volume id_ig_ratio n_at c_at
## 1 CNF/RGO/moOxNy 3.45 0.046 1.45 2.10 18.45
## 2 sulfur-doped graphene foam (SGF) 37.00 0.500 1.28 0.00 85.60
## 3 sulfur-doped graphene foam (SGF) 37.00 0.500 1.28 0.00 85.60
## 4 sulfur-doped graphene foam (SGF) 37.00 0.500 1.28 0.00 85.60
## 5 sulfur-doped graphene foam (SGF) 37.00 0.500 1.28 0.00 85.60
## 6 moS2-moO2/G 3.45 0.170 1.22 1.88 30.75
## o_at electrolyte_formula electrolyte_cond electrolyte_conc
## 1 44.31 H2SO4 7 1
## 2 9.10 KOH 6 6
## 3 9.10 KOH 6 6
## 4 9.10 KOH 6 6
## 5 9.10 KOH 6 6
## 6 52.48 KOH 6 2
## cell_conf rct_na rs_na n_at_na
## 1 three-electrode system TRUE FALSE FALSE
## 2 two-electrode system FALSE FALSE FALSE
## 3 two-electrode system FALSE FALSE FALSE
## 4 two-electrode system FALSE FALSE FALSE
## 5 two-electrode system FALSE FALSE FALSE
## 6 three-electrode system TRUE FALSE TRUE
nrow(df)
## [1] 925
Poniżej wyświetlono kilka metryk dla każdej kolumny ze zbioru
summary(df)
## pw_low pw_up pw current_density
## Min. :-1.1000 Min. :-0.2000 Min. :0.4000 Min. : 0.05
## 1st Qu.:-0.3000 1st Qu.: 0.4000 1st Qu.:0.6000 1st Qu.: 1.00
## Median : 0.0000 Median : 0.6000 Median :0.8000 Median : 2.00
## Mean :-0.2335 Mean : 0.6302 Mean :0.8628 Mean : 5.80
## 3rd Qu.: 0.0000 3rd Qu.: 0.8000 3rd Qu.:1.0000 3rd Qu.: 5.00
## Max. : 0.2000 Max. : 3.5000 Max. :3.5000 Max. :200.00
## capacitance ssa rct rs
## Min. : 1.4 Min. : 8.896 Min. : 0.080 Min. : 0.200
## 1st Qu.: 150.0 1st Qu.: 105.000 1st Qu.: 1.100 1st Qu.: 0.500
## Median : 262.5 Median : 159.970 Median : 2.010 Median : 0.760
## Mean : 419.3 Mean : 427.082 Mean : 2.761 Mean : 1.351
## 3rd Qu.: 509.9 3rd Qu.: 577.100 3rd Qu.: 3.520 3rd Qu.: 2.000
## Max. :3344.1 Max. :2400.000 Max. :24.200 Max. :17.500
## electrode_conf pore_size pore_volume id_ig_ratio
## Length:925 Min. : 0.53 Min. :0.0200 Min. :0.12
## Class :character 1st Qu.: 6.20 1st Qu.:0.0460 1st Qu.:0.94
## Mode :character Median : 6.20 Median :0.2160 Median :1.09
## Mean :10.29 Mean :0.3702 Mean :1.09
## 3rd Qu.:10.80 3rd Qu.:0.4610 3rd Qu.:1.18
## Max. :44.13 Max. :2.3500 Max. :2.90
## n_at c_at o_at electrolyte_formula
## Min. : 0.00 Min. : 1.40 Min. : 1.90 Length:925
## 1st Qu.: 1.88 1st Qu.:30.75 1st Qu.:12.83 Class :character
## Median : 1.88 Median :54.18 Median :18.50 Mode :character
## Mean : 2.49 Mean :55.97 Mean :22.49
## 3rd Qu.: 2.10 3rd Qu.:81.30 3rd Qu.:26.70
## Max. :23.82 Max. :98.10 Max. :54.28
## electrolyte_cond electrolyte_conc cell_conf rct_na
## Min. :1.000 Min. :0.100 Length:925 Mode :logical
## 1st Qu.:6.000 1st Qu.:1.000 Class :character FALSE:139
## Median :6.000 Median :1.000 Mode :character TRUE :786
## Mean :5.818 Mean :2.698
## 3rd Qu.:7.000 3rd Qu.:6.000
## Max. :8.000 Max. :6.000
## rs_na n_at_na
## Mode :logical Mode :logical
## FALSE:153 FALSE:235
## TRUE :772 TRUE :690
##
##
##
Poniżej umieszczono rozkłady wartości wszystkich wartości
numerycznych
Widzimy, że wartości w zbiorze są szeroko rozproszone, a często pojawia się w nich szum informacyjny, gdzie maksymalnie kilka wartości odbiega od ogólnego trendu. Spróbujmy odfiltrować ten szum o 10% najniższych i najwyższych wartości.
Oprócz kolumn z wartościami numerycznymi, możemy też skategoryzować atrybut “cell_conf” i policzyć jego występowanie. Podobnie można zrobić dla “electrode_conf” i “electrolyte_formula”, lecz liczba ich kategorii jest zbyt duża aby miało sens wyświetlenie ich na wykresach.
Do znalezienia korelacji między zmiennymi wykorzystano macierz korelacji.
Wyszukajmy wszystkie koreleacje najbardziej warte odnotowania. Jako takie uznamy takie pary, których korelacja znajduje się między 0.99 i 0.4 oraz między -0.99 i -0.4.
cor_df_high <- cor_df %>%
filter((abs(correlation) > 0.4 & abs(correlation) < 0.99)) %>%
filter(var1 > var2) %>%
arrange(correlation)
kable(cor_df_high)
| var1 | var2 | correlation |
|---|---|---|
| o_at | c_at | -0.7010226 |
| capacitance | c_at | -0.4553412 |
| pw | o_at | -0.4164487 |
| pw | c_at | 0.4714214 |
| pw_up | pw | 0.5970643 |
| pw_up | pw_low | 0.6370278 |
Analizując znalezione korelacje możemy pominąć korelacje pw_up z pw_low oraz pw_up z pw, gdyż wydają się dosyć naturalne. Jako najciekawsze uznałem zależności między: - capacitance a c_at, - pw a c_at.
Obydwa wykresy pokazują ciekawe właściwości węgla. Możemy z nich wyczytać, że wraz ze zwiększeniem zawartości węgla w materiale jego pojemność właściwa spada. Z drugiej strony rośnie okno w którym materiał pozostaje stabilny. Sprawdźmy jeszcze zależność między pojemnością a oknem stabilności, aby sprawdzić czy one też nie są ze sobą skorelowane.
Jak widać na wykresie, korelacja między nimi również występuje. Natomiast nakreślony trend dla powyżej 2000 F/g jest bardzo wątpliwy z uwagi na stałą wartość 0.5V dla tych pojemności, która może wynikać z metody uzupełniania pustych wartości lub samego zbioru danych.
Przyjrzyjmy się jeszcze pojemności w relacji do składu atomów.
Z wykresu wyżej widzimy, że większy udział tlenu wspiera do pewnego stopnia pojemność materiału. Tak jak wcześniej zauważyliśmy, wpływ węgla jest negatywny na pojemność, a jego odbicie na końcu prawdopodobnie jest jedynie błędem z powodu małej próbki. Ostatnim wnioskiem jest bardzo niski udział azotu we wszystkich materiałach, a jego zmiana ma marginalny wpływ na pojemność.
W tej części spróbujemy przewidzieć pojemność materiału na podstawie tego zbioru danych.Na początku spróbujmy znaleźć wartościowe parametry. Skorzystamy z poprzedniego kodu. Przyjmijmy, że kolumny warte zachowania to są te których korelacja jest powyżej 0.1 lub poniżej -0.1 (z pominięciem pw_low i pw_up, gdyż zakładam ich dużą korelację z pw).
Spróbujmy zamodelować pojemność regresją liniową.
cor_df_cap <- cor_df %>%
filter(var1 == 'capacitance' & var2 != 'capacitance' & abs(correlation) > 0.1) %>%
arrange(correlation)
kable(cor_df_cap)
| var1 | var2 | correlation |
|---|---|---|
| capacitance | c_at | -0.4553412 |
| capacitance | pw | -0.3384832 |
| capacitance | ssa | -0.2209595 |
| capacitance | pw_up | -0.1047333 |
| capacitance | pore_volume | -0.1010717 |
| capacitance | id_ig_ratio | 0.1026857 |
| capacitance | pw_low | 0.1752674 |
| capacitance | o_at | 0.3204094 |
df_capacitance <- df %>%
select(capacitance, c_at, pw, ssa, pore_volume, id_ig_ratio, o_at)
inTraining <-
createDataPartition(
y = df_capacitance$capacitance,
p = .75,
list = FALSE)
training <- df_capacitance[ inTraining,]
testing <- df_capacitance[-inTraining,]
ctrl <- trainControl(
method = "repeatedcv",
number = 2,
repeats = 5
)
fit <- train(capacitance ~ .,
data = training,
method = "lm",
trControl = ctrl
)
rfClasses <- predict(fit, newdata = testing)
metrics <- postResample(pred = rfClasses, obs = testing$capacitance)
metrics
## RMSE Rsquared MAE
## 387.8635335 0.2515056 268.3858543
results_df <- data.frame(
Actual = testing$capacitance,
Predicted = rfClasses
)
min_val <- min(c(results_df$Actual, results_df$Predicted))
max_val <- max(c(results_df$Actual, results_df$Predicted))
p <- ggplot(results_df, aes(x = Actual, y = Predicted)) +
geom_point(alpha = 0.6, color = "#1f77b4") +
geom_abline(intercept = 0, slope = 1, linetype = "dashed", color = "red", linewidth = 1) +
geom_smooth(method = "lm", color = "darkgreen", se = FALSE) +
coord_fixed(xlim = c(min_val, max_val), ylim = c(min_val, max_val)) +
labs(
title = "Prognozowana vs. rzeczywista pojemność - regresja liniowa",
x = "Rzeczywista pojemność (F/g)",
y = "Przewidywana pojemność (F/g)"
) +
theme_minimal()
ggplotly(p)
Na wykresie widzimy 2 linie - zielona to nasz wytrenowany model, a czerwona to idealna regresja dla tego zbioru danych. Jak można na wykresie zauważyć, nasz model kompletnie nie radzi sobie z większymi pojemnościami. Spróbujmy skorzystać z algorytmu Random Forest Regression.
head(training)
## capacitance c_at pw ssa pore_volume id_ig_ratio o_at
## 2 367 85.60 1.0 537.0 0.50 1.28 9.10
## 3 338 85.60 1.0 537.0 0.50 1.28 9.10
## 4 283 85.60 1.0 537.0 0.50 1.28 9.10
## 6 872 30.75 0.5 168.2 0.17 1.22 52.48
## 7 143 54.18 0.6 385.0 0.10 1.00 12.83
## 10 483 54.18 0.6 385.0 0.31 1.21 12.83
rf_fit <- train(capacitance ~ .,
data = training,
method = "rf",
trControl = ctrl,
tuneLength = 3
)
rf_predicts <- predict(rf_fit, newdata = testing)
rf_metrics <- postResample(pred = rf_predicts, obs = testing$capacitance)
results_df <- data.frame(
Actual = testing$capacitance,
Predicted = rf_predicts
)
min_val <- min(c(results_df$Actual, results_df$Predicted))
max_val <- max(c(results_df$Actual, results_df$Predicted))
p <- ggplot(results_df, aes(x = Actual, y = Predicted)) +
geom_point(alpha = 0.6, color = "#1f77b4") +
geom_abline(intercept = 0, slope = 1, linetype = "solid", color = "red", linewidth = 1) +
geom_smooth(method = "lm", color = "darkgreen", se = FALSE, linewidth = 0.8) +
coord_fixed(xlim = c(min_val, max_val), ylim = c(min_val, max_val)) +
labs(
title = "Prognozowana vs. rzeczywista pojemność - RF",
x = "Rzeczywista pojemność (F/g)",
y = "Przewidywana pojemność (F/g)"
) +
theme_minimal()
ggplotly(p)
Jak widzimy na wykresie, algorytm Random Forest poradził sobie o rząd wielkości lepiej aniżeli regresja liniowa. Spójrzmy na sam koniec jak wyglądają przykładowe wartości Shapleya dla jednego przykładu.
testing[1,]
## capacitance c_at pw ssa pore_volume id_ig_ratio o_at
## 1 680 18.45 0.8 186.3 0.046 1.45 44.31
ive_rf <- individual_variable_effect(rf_fit, data = testing,
new_observation = testing[1,])
plot(ive_rf)